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Using coupled equations for the bosonic and fermionic order parameters, we construct families of 
gap solitons (GSs) in a nearly one-dimensional Bose-Fermi mixture trapped in a periodic optical- 
lattice (OL) potential, the boson and fermion components being in the states of the BEC and BCS 
superfluid, respectively. Fundamental GSs are compact states trapped, essentially, in a single cell of 
the lattice. Full families of such solutions are constructed in the first two bandgaps of the OL-induced 
spectrum, by means of variational and numerical methods, which are found to be in good agreement. 
The families include both intra-gap and inter-gap solitons, with the chemical potentials of the boson 
and fermion components falling in the same or different bandgaps, respectively. Nonfundamental 
states, extended over several lattice cells, are constructed too. The GSs are stable against strong 
perturbations. 



PACS numbers: 03.75.Ss,03.75.Lm,05.45.Yv 



I. INTRODUCTION 

The experimental realization of Bosc-Einstein conden- 
sation (BEC) [1 was followed by the creation of another 
ultra-cold low-density quantum medium, viz., the degen- 
erate Fermi gas (DFG) [2| • Fermi gases cannot be directly 
produced by means of the evaporative cooling, due to 
the Pauli blockade of interactions among spin-polarized 
fermions. However, it has been possible to achieve the 
DFG state using the sympathetic-cooling technique, in 
the presence of another boson or fermion component in 
the gas. The use of this technique has resulted in the ob- 
servation ^M, IH 5 1 and related experimental and 
theoretical % M^UjL 03 studies of trapped DFG-BEC 
mixtures composed of the following species: 7 Li- 6 Li [J], 
23 Na- 6 Li pi, and 87 Rb- 40 K [1,0]. New remarkable exper- 
iments have made use of the Feshbach-resonance tool to 
induce attraction between fermions in spin-up and spin- 
down states, thus creating the Bardeen-Cooper-Schricffcr 
(BCS) superfluid [3 in Fermi gases of 6 Li Q and 40 K 
[l5j atoms. Actually, the attraction in these experiments 
was varied from weak, appropriate to the BCS state per 
sc, to strong, appropriate to the formation of the BEC of 
diatomic molecules, making it possible to explicitly study 
the BCS-BEC crossover, as theoretically predicted in ear- 
lier works and analyzed in more recent ones fl7f fa 
review of experimental results can be found in Ref . IISII , 
and a brief overview of the theory was given in Ref. [l9(). 

Among most fundamental matter-wave patterns ob- 
served in BEC are bright solitons, see original works 
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[20l . |21| and reviews [22| ■ Further analysis has predicted 
the formation of bright solitons in degenerate Bose-Fermi 
[23l . |24| and Fermi- Fermi mixtures [25| , in the presence of 
attractive inter-species interactions, which may overcome 
the effective Pauli repulsion between identical fermions 
(and possible repulsion in the BEC component). 

On the other hand, gap solitons (GSs) have been pre- 
dicted [Hi and further studied [H, 0, [H [H H| in 
BEC with the repulsion between atoms, trapped in the 
periodic potential induced by an optical lattice (OL). GSs 
are localized objects supported by the balance between 
the self-repulsion and negative effective mass provided 
by the OL-induced matter-wave spectrum. The chemical 
potential of the GS must fall in a spectrum's bandgap, 
where delocalized Bloch-wave states cannot exist. 

The prediction of GSs in the BEC Hi| was followed 
by the creation of a GS formed by ~ 250 atoms of 87 Rb 
in a "cigar-shaped" t rap combined with the longitudi- 
nal OL potential [HI, In the experiment, the BEC 
was pushed in a state with the negative effective mass 
by means of acceleration. Another possibility may be to 
add a strong parabolic trap to the OL potential, thus 
squeezing the system into a small region and then slowly 
relaxing the extra trap, to give (a part of) the condensate 
a chance to remain in a squeezed GS state [35| . Nonfun- 
damental states, in the form of a stretched "gap-wave" 
pattern, exten ding over many cells of the OL, have also 
been observed [36[ and interpreted theoretically as a seg- 
ment of a nonlinear Bloch wave confined by sharp fronts 
[37l ] (in fact, similar states were also predicted in the self- 
attractive BEC mi). 

Multicomponent GSs were predicted too. A straight- 
forward possibility is to consider a symmetric model for 
a mixture of two boson species with equal masses and 
coefficients of the coupling to the OL potential, which 
may represent two different hypcrfinc states of the same 
atom, with repulsion between the species. As shown in 
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Ref. [30j . this model gives rise to stable 2D and ID 
gap solitons - first of all, in the most fundamental case 
when the intra-specics repulsion is absent (as may be ad- 
justed by means of the Feshbach-resonance technique), 
and the GSs are supported, as symbiotic [H] states, solely 
by the inter-species interaction. These two-component 
GSs may be of intra-gap and inter-gap types, with the 
chemical potentials of the components falling into a sin- 
gle bandgap, or belonging to different gaps, respectively 
(inter-gap solitons are essentially less stable [301 )■ Also 
studied were three-component GSs in the ID model of a 
spinor BEC trapped in the OL [32l |. 

In this work we aim to study the possibility of the 
formation of fundamental GSs in a one-dimensional (ID) 
model of the intrinsically repulsive superfluid Bose- Fermi 
mixture (SBFM), whose boson and fermion components 
are assumed to be in the BEC and BCS-superfluid 
phases, respectively. Fundamental GSs, unless taken very 
close to edges of the corresponding bandgap, feature a 
compact shape, being essentially localized in a single well 
of the OL potential, without tangible undulating tails. 
Obviously, a localized pattern cannot be a ground state of 
the gas with repulsive interactions. Nevertheless, the GSs 
in intrinsically self-repulsive BEC have been found to be 
stable against small or moderate perturbations (2?| (i.e., 
they are, generally speaking, metastable states), and, as 
mentioned above, they have been created in experiment 
[33l . HH. We also show that, in addition to the compact 
fundamental GSs, one can find solutions for nonfunda- 
mental solitons, which occupy several OL sites and are 
dynamically stable too. 

It is well known that a dilute BEC is very accurately 
described by Gross-Pitaevskii equation (GPE) for the 
mean- field wave function, which, in particular, provides 
for a reliable model for the dynamics of dark and bright 
solitons, vortices, onset of collapse, etc. 0, [22J]. In 
the Bose-Fermi mixture, the bosonic component is also 
treated by means of the GPE (with an extra term tak- 
ing into regard collisions with fermions), whereas the 
fermionic part should be, in principle, described by a set 
of Schrodinger equations for individual atoms (4(| ■ How- 
ever, with the increase of the number of fermions such ab 
initio treatment becomes much too complicated to im- 
plement. On the other hand, a simplified model of the 
BCS superfluid of the Ginzburg-Landau type is widely 
accepted [Hj]. This approach describes the fermion su- 
perfluid by a single complex order parameter, which may 
be interpreted, at a semi-microscopic level, as the wave 
function of composite bosons (Cooper pairs). The model 
based on nonlinearly coupled equations for the boson and 
fermion order parameters has been used to study various 
phenomena in Bose-Fermi mixtures, such as vortex for- 
mation, mixing-demixing transitions, and the onset of 
collapse [ll|, [4l[ [42|] . We will use these equations for the 
prediction of compact GSs in the SBFM, which seems 
quite a relevant approach to the description of these 
simple-shaped objects. 

The so derived model for the SBFM in three dimen- 



sions amounts to a system of nonlinear Schrodinger equa- 
tions, with the cubic nonlinearity for the bosons, and 
nonlinear term of power 7/3 for the fermions (43|, |44| ; 
the two equations are coupled by cubic terms. In the 
presence of a strong transverse trap, this system can be 
reduced to coupled ID equations. We look for GS solu- 
tions and analyze their stability within the framework of 
the ID system. The simplicity of the Gaussian-like pro- 
files of the GSs suggests to apply the variational approx- 
imation (VA) [IH to them, which produces very accurate 
results, if compared to numerical findings. 

It is relevant to mention that states of the GS type 
in a Bose-Fermi mixture were also investigated in Ref. 
[io| , using a microscopic description of the fermion com- 
ponent. Unlike the compact GSs studied in the present 
work, the solitons found in Ref. [I(| typically occupy 
several OL cells, and feature a sophisticated shape with 
undulating tails. Another difference is that Ref. [4(| 
was chiefly dealing with the attractive Bose-Fermi inter- 
action, whereas in our system both the Bose-Fermi and 
Bose-Bosc interactions could be repulsive as well as at- 
tractive. 

The paper is organized as follows. In Sec. II we de- 
rive the basic equations from the underlying Lagrangian 
density. Section III is dealing with the VA based on the 
Gaussian ansatz for the wave functions. In Sec. IV we 
report numerical results for intra-gap and inter-gap GS 
families in the first two bandgaps, and compare them 
to the VA predictions, concluding that the agreement 
between the numerical and variational methods is very 
good. We also verify the stability of the GSs by means 
of direct simulations. Section V summarizes the work. 



II. COUPLED EQUATIONS FOR THE 
SUPERFLUID BOSE-FERMI MIXTURE (SBFM) 

A. Three-dimensional equations 

We consider a degenerate mixture of Nb bosons and 
Np fermions with spin ^ , whose masses are tub and m p , 
at zero temperature. The boson and fermion compo- 
nents of the mixture are assumed to be condensed into 
the BEC and a BCS superfluid. respectively Accord- 
ingly, a weakly attractive underlying (hidden) interac- 
tion between spin-up and spin-down fermions is implied, 
while the explicit fermion-boson and boson-boson inter- 
actions are repulsive, which is the most natural case for 
any atomic gas. 

The derivation of the model starts with the effective 
Lagrangian density for the SBFM. 

C = Cb + Cf + Cbf, (1) 
with the usual expression for the boson component [lj, 

ih ( dt[>* B di> B \ h 2 2 
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where iftg is the mean- held order parameter (single-atom 
wave function) of the BEC, the density of the boson 
atoms being iig = \iPb\ ■> and Vs(r) is the trapping 
potential for the bosons. The strength of the boson- 
boson repulsive interaction is Gb = Awh^aBB/mB, with 
cibb > the respective scattering length. 

To derive the Lagrangian density of the fermion com- 
ponent, we use the known energy density of the BCS 
superfluid (SF) [H, S3, SI , 

£sf = (3/5)n F e F , (3) 

where £p = h 2 kp/(2m) is the Fermi energy, hkp is the 
Fermi momentum, and np the density of the fermion 
atoms. This energy density was first derived by Lee and 
Yang [47} in the weak-coupling BCS limit. Modifications 
to this expression for the description of the BCS-BEC 
crossover (with the gradual increase of the strength of the 
attraction between spin-up and spin-down fcrmions) have 
also been considered [48[. With regard to the pairing of 
the fcrmions with opposite spin orientations, the total 
fermion density is given by np = 2(2-7r)~ 3 ] J Q F 47rfc 2 dfc = 

4 /2 M 3/2 , with A = (3TT 2 ) 2 / 3 h 2 /(2m F ), which yields the 
expression for the Fermi energy in terms of the density: 
ep = A(uf /2) 2 / 3 . Then, the energy density in Eq. © 

is expressed as £sf = {3/5)An 5 J 3 , and the respective 
fermionic Lagrangian density in Eq. |T]) becomes 



ih ( dip* F 
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where ipp is the above-mentioned complex order param- 
eter of the fermionic BCS condensate, which determines 
its density, np — \iJjf\ 2 , VfC 1 ") is the potential trapping 
the fermion atoms, and rhp is an effective mass of the 
superfluid flow; (it is natural to expect that in a weakly 
coupled superfluid rhp is close to 2mp.) 

Lastly, the Bose-Fermi repulsive interaction is ac- 
counted for by the corresponding term in the Lagrangian 
density, 
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where Gbf — 2irh 2 a,BF /ran, with the respective reduced 
mass, nriR = m,BVHp j [tub + nip), and scattering length 
cibf- Substituting expressions ©, ©, and © in full 
Lagrangian density ([1]), we derive the following Eulcr- 
Lagrange equations, cf. Ref. [l2l ]: 



^ d h 2 

-in— 

dt 2ms 



VI + Vb{t)+G b \M' 



+G bf \^f\ 2 ] B (r,t) = O, 

d n 2 v 2 



dt 2m p 



V F (v)+A\ip t 
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+G BF \iPb\ 2 ^jr(r,i)=0, 
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(7) 



which are supplemented with normalizations 



\i>B,F(x,V,z)\ dxdydz = N 



B.F 



(8) 



[these two norms are dynamical invariants of Eqs. ^ 
and ©]. 

Dynamical equations ^ and ([7]) can also be de- 
rived in a different physical context, viz., in a DFG 
in the hydrodynamic limit [49| . Indeed, using the 
known energy density for the DFG, and a kinetic en- 
ergy term h 2 / {2mp)\ V-0f| 2 for the hydrodynamic flow 
in the Lagrangian density (here, niF is the effective 
mass of the hydrodynamic flow in the DFG), one ar- 
rives at Eqs. ([6]) and {2J with the difference that 
the term A = (3-7T 2 ) 2 / 3 h 2 / \2mp) therein is replaced by 
(6n 2 ) 2 ^ 3 h 2 / (2m.F)', the difference in the factor of 2 is due 
to the formation of Cooper pairs in the superfluid. Hence 
the present analysis and conclusions for the stationary 
GSs should also apply to the DFG. However, this anal- 
ogy may be questionable for nonstationary states, as the 
above-mentioned description of the DFG does not define 
the phase of the fermion order parameter. 



B. Reduction to one-dimensional equations 

As said above, the presence of a strongly elongated 
(cigar-shaped) trap suggest to reduce the 3D equations 
to a ID system. In the usual GPE with the cubic non- 
linearity, the reduction of the full 3D equation to its ID 
counterpart was performed, in different ways, in Rcfs. 
[ITj, HH- I n the simplest situation, the reduction of the 
3D equation to one dimension starts with the factoriza- 
tion of the 3D wave function, 



iPb,f(x, y, z, t) = <I>b,f{x, t) exp —i (cjx) 



B.F 



t - 



y 2 + z 2 



2 (&ho) b F 

(9)' 

where the second multiplier represents the ground state 
of the transverse 2D harmonic oscillator, with (ah a ) B F = 



h/ ms,F (w_l) 



B.F 



the harmonic-oscillator lengths for 



the bosons and fermions, and (uj_) B F the corresponding 
transverse trap frequencies. The normalization of func- 
tions 4>b,f(x), which were introduced in Eq. ([9]), must 
comply with the underlying 3D normalization conditions 
©. The substitution of ansatz © in Eqs. © and © 
and averaging of the 3D equation in the transverse plane, 
(y,z), lead to the following ID equations: 
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where we have introduced the periodic OL potential with 
period A/2 and amplitudes —£b,f, acting on the Bose 
and Fermi components (the amplitudes are taken neg- 
ative to set a local minimum of the potential at point 
x = 0, where the center of the solitons will be placed). 
Below, results are presented assuming equal boson and 
fermion masses and equal constants of the coupling to 
the OL: m F = n^B = tn, cb = £f = £■ For a 
Bose-Fermi superfluid mixture this situation may approx- 
imately correspond to the existing 40 K- 87 Kb mixture Q, 
where m F ~ 2m p 

By means of the rescalings (here we set m = tuf = "niB 
and €b = e_F = e), 
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Eqs. (JTQ1) and (fTTj) are cast in the following dimcnsionless 
form, 

,9<t>B 190 S 2 

+ .9b|<Pb| 0b 



dt 2 dx 2 

+ 3sf|0f| 2 0b - ecos(2x)0 s , (13) 



190F 



+ .g F |<M 4/3 <^ 



dt 2 ,9a; 2 

+ .9bb|</>b| 2 <?V - ecos(2a;)0 F , (14) 

with .g B = a BB A/(7ra 1 2 io ), g B F = a B F^/(^al Q ) and 
g F = (3/10)(3A 2 /47rog o ) 2/3 (tildes are dropped here), the 
rescalcd ID wave functions being subject to normaliza- 
tion conditions 



\0b,f(x)\ 2 dx = N. 



B.F- 



(15) 



Stationary solutions to Eqs. (|T3|) and (|T4|) arc looked 
for in the ordinary form, 4>b,f = cxp (— i/iB.Ft) ub,f(x), 
with real chemical potentials \ib and /if, and real func- 
tions ub.f(x) obeying equations (with u" = d 2 u/dx 2 ) 



u B u B 



flpUp 



gBU 3 B + gBFUpus - ecos (2x) U B , 



— +9FUp 
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(16) 

9bfU b Uf — e cos (2x) up. 

(17) 



Our main objective is to construct families of fundamen- 
tal GS solutions to Eqs. (fTo) and (fl7|) . with yU_B, F belong- 
ing to the first two finite bandgaps in the model's linear 
spectrum [the bandgaps are induced by the dimcnsionless 
potential, —ecos (2a;), common to both components]. 



III. VARIATIONAL APPROXIMATION (VA) 

To derive the VA [45[ for the fundamental GS solutions 
to stationary equations (fT6|) and (TiT]) . we use the obvi- 
ous fact that these equations, along with normalization 
conditions (|15[) . can be derived from Lagrangian 
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Hbu b 



UpUp 



1 3 

+ ecos(2x)(u% + up) - -gBU% - -gpu 1 ^ 3 
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where [i B and up play the role of Lagrange multipliers 
associated with conditions (fTS")) . Aiming to find solitons 
with a compact profile, we adopt the Gaussian ansatz 

m, 



ub,f{x) 



1/4 INb,f^b,f 
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■ cxp 



2w* B<F 
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where the variational parameters are soliton's widths 
wb.f arid reduced norms ^b,f, as well as chemical poten- 
tials hb,f- The substitution of the ansatz in Lagrangian 
(p|) yields 



L = hbN b K b + u f NfKf 
+ eN B K B e- w z + eN F K F e 
9f 
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The first pair of the variational equations following from 
Eq. (HOI), dL/d/iB.N = 0, yield K B = K F = 1 [which im- 
plies that the norms of the two components of ansatz 
(fill)) arc Nb and N F , in accordance with Eq. (fl~5"|) ]. 
Therefore, K F = N F = lis substituted in other equa- 
tions below. The other pair of the variational equations, 
dL/du'B,N = 0, can be written as 
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(22) 

The remaining equations, dL / 9&b,f = 0, yield expres- 
sions for chemical potentials /is, F of the boson and 
fermion components of the fundamental GS as functions 
of the interaction coefficients, g Bl gF, and g B F- 
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(23) 
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(24) 

In the limit of g B>F , g B F -> 0, Eqs. |21]) and ([22]) 
reduce to a single equation, 

Aew% >F exp (-«;|,,f) = 1, (25) 

which has real solutions for e > e 2 /8 f=a 0.92. Formally, 
solutions of Eq. (f2"5")) predict localized solutions to Eqs. 
(fl"6|) and (| 1 7[) in the linear system, while, as is com- 
monly known, the linear Schrodinger equation with a 
periodic potential does not support any localized solu- 
tion. However, the actual meaning of the limit case of 
9b,f, 9bf — > is that, with solutions to Eq. (|25|) substi- 
tuted in Eqs. (|2U|) and (|2"4")h the latter equations predict 
a cutoff value, (fio) 1 , for both chemical potentials, /ib 
and hf, below which fundamental GSs do not exist, i.e., 
as a matter of fact, (/io)i is the left edge of the first (low- 
est) bandgap. A numerically computed position of the 
left edge agrees very well with the variational prediction. 
For instance, for e = 5 (numerical results are displayed 
below for this value of e), the VA yields (a*o) 1 ~ —2.894, 
which almost exactly coincides with its numerically found 
counterpart, (Ato)i num ' ~ —2.893. 

IV. NUMERICAL RESULTS 

A. Families of fundamental gap solitons 

For the numerical solution, Eqs. (|13p and Q14D were 
discretized using the Crank-Nicholson scheme [52| and 
solved (in real time) until the solution would converge 
to a stationary wave form that provides a real solution 
to Eqs. (flB|) and (fTT|) . The numerical solutions was car- 
ried out with time step 0.0005 and space step 0.025, in 
domain —20 < x < 20. To find the GSs, the simula- 
tions started with an initial configuration chosen as the 

ground state, 4>(x) = (y/2c/ir) 1 / i exp (^—^c/2x 2 ^j, of the 

quantum-mechanical linear oscillator with potential cx 2 
(with large c, typically 5 to 200). The OL potential was 
slowly introduced in the course of numerical simulation, 
while the strong harmonic potential kept the localized 
state squeezed into a single cell of the OL. After reaching 
a stationary bound state in the nonlinear equation con- 
taining the combined OL and harmonic potential, the 
latter term was slowly switched off, allowing the estab- 
lishment of a well-defined stationary GS (a similar ap- 
proach to the creation of GSs in the single-component 
BEC was proposed in Ref. [35|). 

Families of the fundamental GSs predicted by the VA 
may be characterized by dependences <7-B,-F (/■*)• The same 
dependences were found from the numerical solution out- 
lined above. In this subsection, we display numerical 
results obtained for compact GSs, without conspicuous 
"tails" . 



The effective nonlinearities, gsNs and gpN F , for the 
boson and fermion components of the GS families are 
shown, against the corresponding chemical potentials, in 
Figs. [IJa) and (b), fixing e = 5 and gBF = 0.004. In 
panel (a), we also fix gp = 0.05 and Nb = 1000, while 
the bosonic-nonlincarity coefficient, gB, varies from to 
0.015, and the results are displayed for different fixed 
values of Np (in particular, Np = corresponds to the 
ordinary GSs in the self-repulsive BEC). The families of 
variational solutions displayed in Fig. Ufa) commence at 
/i w —2.894 (as said above, it almost exactly coincides 
with the actual left edge of the first bandgap). In panel 
(b) of Fig. [Q we fix g B = 0.01 and Np = 0.05, while the 
fermion-nonlinearity coefficient, g F , varies from to 0.1, 
the results being shown for several different fixed values 
of Nb (the case of Nb = corresponds to pure fermionic 
GSs in the BCS superfluid). The agreement between the 
VA and numerical results is very good, with a caveat that 
the VA ignores the presence of the Bloch band separating 
the two gaps (i.e., the VA-generated solution branches 
extend across the band). 

In Fig. [21 typical shapes of the numerically found 
GSs are compared with their counterparts predicted by 
the VA at the following sets of parameter values: (a) 
Np = 1000, N B = 500, g B = 0.01; (b) N F = 200, N B = 
1000,3b = 0.001; (c) N F = 100, N B = 700, g B = 0.01; 
(d) N F = 100, A B = 500, g B = 0.01. In all the four 
cases, gp = 0.05 and gBF = 0.004 are fixed. As written 
in the caption to Fig. [2j these parameters are chosen to 
represent four different cases, as concerns the location of 
the two chemical potentials, (is and MFj m the two low- 
est bandgaps of the OL-induccd linear spectrum: panels 
(a) and (d) represent the two-component solitons of the 
intra-gap type, while (b) and (c) represent inter-gap soli- 
tons, in terms of Ref. [3Q|. The agreement between the 
numerical and variational results is, again, quite note- 
worthy. 

In Fig. [2] all the nonlinearities are taken as positive 
(repulsive) . Another problem of interest is to investigate 
if GSs in the first two bandgaps can be found in both 
components when one of the nonlinearities (gB or gBF) 
becomes attractive. (We recall that when gsF is suf- 
ficiently attractive, bright Bose and Fermi solitons can 
be generated without the OL potential [24|. Such bright 
solitons will not have the chemical potential in the first 
two bandgaps and hence are not GSs.) We find that GSs 
in the first two bandgaps can indeed be found when ei- 
ther gs or gBF turns negative. We show typical results 
for this case in Fig. [3] (a) and (b), where we illustrate 
that by adjusting the parameters such GSs can be placed 
in the first or the second bandgap. These are inter-gap 
solitons. In addition intra-gap solitons can also be cre- 
ated in this case (not illustrated in Fig. [3]) The solitons 
shown in Fig. [3] are also stable, despite the fact that one 
nonlincarity coefficient changed its sign to an attractive 
value. 
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FIG. 1: (Color online) The effective nonlinearity coefficients, 
gsNs for bosons (a) and gFNp 3 for fermions (b), versus the 
chemical potentials of the boson and fermion components, fiB 
and fip, respectively, for families of the fundamental compact 
gap solitons in the first and second bandgaps of periodic po- 
tential ecos(2a-) (numerical results are presented in this and 
other figures for e = 5) . As indicated in the figure, continuous 
curves and chains of symbols display, respectively, results ob- 
tained from the numerical solution ("num") and variational 
approximation ( "var" ) . Shaded vertical areas represent Bloch 
bands which separate the gaps. 



B. Nonfundamental solitons 




FIG. 2: (Color online) Profiles of the probability density, 
\4>b,f(x)\ 2 [subject to normalization \(f>B,F(x)\ 2 = 1] in 
the boson ( "bos" ) and fermion ( "fer" ) components of the fun- 
damental gap solitons, in four different cases with (a) both 
components of the soliton having their chemical potentials, 
\ib and fiF, in the second bandgap, (b) ub and fiF in the 
first and second gaps, respectively; (c) vice versa to (b); (d) 
both components the first bandgap. In terms of Ref. [30], 
the two-component gap solitons in panels (a) and (d) are of 
the intra-gap type, while the ones in (b) and (c) are inter- 
gap solitons. The thin dotted lines, here and in other figures, 
represent the OL potential, oc — cos(2a:). 



The fundamental GSs presented by Figs. Q] and [5] are 
compact objects with a nearly-Gaussian shape, trapped 
in a single cell of the OL. On the other hand, it is well 
known [26[ that the GSs in self-repulsive repulsive BEC 
frequently demonstrate structures extended over many 
cells. Stable nonfundamental solitons featuring similar 
profiles can be found in the present model too, see ex- 
amples in Fig. 3J They may be considered as bound 
complexes built of a central GS and additional pulses lo- 
cated in adjacent cells of the OL. Similar to the patterns 
shown Fig. [4] which extend over three lattice cells, at 
smaller values of the chemical potentials it is possible 
to find dynamically stable localized structures which are 
stretched still wider. 



C. Stability 

The way the fundamental and nonfundamental GSs 
have been found, as a limit form to which solutions of un- 
derlying dynamical equations fT3|) and ([14")) relax at large 
values of t, strongly suggests that all these solutions are 
stable. Their stability was additionally tested by adding 
strong perturbations to them. It was thus verified that 
the entire families of the GSs presented above are sta- 
ble, as well as the nonfundamental solitons displayed in 
Fig. |U For example, the nonfundamental soliton labeled 
III in Fig. 21 after it was obtained as a steady state to 
which the evolution of a certain initial configuration con- 



7 



1.2 
0.8 
£ 0.4 



N 



N, 



400 
i B = 1000 
g B = -0.001 
g F = 0.05 
9bf = 0.004 
U B = -2.72 
(?=1.79 




bos (var) 
fer (var) 
bos (num) 
fer (num) 



(a) -3 -2 



1.2 



0.8 



£ 0.4 



N 



N, 



100 
i B = 700 
Or = 0.01 
g F = 0.05 
ge F = -0.002 
(l B = 1 .33 
(l F = -2.92 




bos (var) 
fer (var) 
bos (num) 
fer (num) 



(b) 



-2 



FIG. 3: (Color online) Same as in Fig. [2]with (a) intraspecies 
attraction among bosonic atoms corresponding to a nega- 
tive gB and (b) interspecies attraction between bosons and 
fermions corresponding to a negative Qbf- In (a) the bosonic 
GS is in the first bandgap and fermionic GS is in the second 
bandgap, in (b) vice versa to (a). 



verged, was suddenly multiplied by a perturbing factor, 
4>b,f{x) — > 1.14>b.f{%)- The perturbed solution demon- 
strates oscillations without any onset of instability, as 
shown in Fig. [5] 

A typical example of the stability test for fundamental 
GSs is presented in Fig. II V C| for the soliton shown in 
Fig. [Ha). At t = 10, it was also "shaken", multiplying 
it by factor 1.1. The subsequent evolution makes the 
absence of any instability evident. 



V. CONCLUSION 

We have studied the possibility of generation of one- 
dimensional GSs (gap solitons) in a Bose-Fermi mixture 
trapped in a periodic OL potential, assuming that the 
bosons are in the BEC phase, and fermions are condensed 
into the BCS superfluid. A dynamical model of the mix- 
ture was derived from the system's Lagrangian. The re- 
sultant coupled equations feature the cubic nonlinearity 
in the boson component, and the nonlinear term of power 
7/3 in the equation for the fermion order parameter, with 
a cubic coupling between the equations, all the nonlin- 
caritics being repulsive. 

Families of solutions for compact GSs, trapped, essen- 
tially, in the single lattice cell, were obtained by means of 
the VA (variational approximation) based on the Gaus- 
sian ansatz, and in the numerical form, in the two lowest 
bandgaps of the OL-induced spectrum. For a wide range 
of values of the interaction strengths and numbers of bo- 
son and fermions in the mixture, the VA predictions for 
the GS profiles and their chemical potentials agree very 
well with the respective numerical results. The GS fam- 
ilies include both intra- and inter-gap solitons, with the 
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FIG. 4: (Color online) Profiles labeled I, II, III [normalized 
according to \4>b,f{x)\ 2 — 1] represent the bosonic (a) and 
fermionic (b) components of three distinct stable nonfunda- 
mental gap solitons, extended over several lattice cells, for 
the same parameters (except for [ib,f) as in Fig. [2ja). The 
corresponding values of the chemical potentials are indicated 
in the panels. 




FIG. 5: (Color online) The evolution of the boson (a) and 
fermion (b) components of the nonfundamental gap soliton, 
labeled III in Figs. g] with fi B = -0.80 and fi F = -0.60. 
A rather strong perturbation was applied to the soliton at 
t = 10, in the form of 4>b,f{x) — > \.14>b,f(x). The perturbed 
soliton remains stable as long as the simulation was running. 




FIG. 6: (Color online) The same stability test as shown in 
Fig. [5] for the fundamental soliton from Fig. [2] 
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chemical potentials of the boson and fcrmion components 
falling, respectively, in the same or different bandgaps. 
In addition to the fundamental GSs, we have also found 
(in the numerical form) nonfundamental localized struc- 
tures, extended over several lattice sites. The entire fami- 
lies of the fundamental GSs, as well as the nonfundamen- 
tal extended states, were checked to be robust against 
strong perturbations. 

This work can be naturally extended in various direc- 
tions. One possibility is to look for GSs in 2D and 3D 
models of the same type, study the possibility of the exis- 
tence and stability of the two-component boson-fermion 
gap solitons and investigate their properties. Another 
perspective direction would be to study gap solitons in 
a boson-fermion mixture with mixed repulsive-attractive 
interactions in a systematic form, while in this paper we 
gave only a few examples. One can also seek for mixed 
bound states, in which the boson component belongs to 
the semi-inhnite gap, while the fermionic chemical po- 
tential still belongs to a finite bandgap and vice versa. 
Lastly, it may happen that the coefficients which couple 
the bosons and fermions to the OL have opposite signs 



(one species is red-detuned, while the other one feels blue 
detuning) , that may give rise to nontrivial stability limits 
for two-components GSs. 

Here we used a set of mean- field equations for the study 
of GSs in Bose-Fermi mixtures. A rigorous treatment 
of this problem might be performed using a many-body 
Slater-determinant wave function (40| . as. for example, 
in the case of many-electron scattering [53j . 
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